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1»  INTRODUCTION.  In  its  application  to  the  solutions  of  engineering 
problems,  the  finite  element  discretization  has  been  implemented  almost, 
exclusively  to  the  spatial  dimensions.  For  dynamic  or  time-dependent 
problems  whose  solutions  as  functions  of  time  are  of  interest,  a step- 
by-step  procedure  of  finite  difference,  i.e.,  the  quasi-static  approach 
is  usually  employed.  The  answer  to  the  question  why  the  time  dimension 
has  not  been  treated  equally  with  the  spatial  variables  in  the  finite 
element  discretization  must  be  related,  in  part  at  least,  to  the 
development  of  variational  methods,  since  the  finite  element  procedure 
can  be  viewed  most  readily  as  an  extremizing  sequence  associated  with  a 
variational  statement.  While  there  are  numerous  variational  principles 
for  boundary  value  problems,  few  exist  for  initial  value  problems.  Like 
many  problems  involving  nonconservative  forces,  the  difficulty  appears  to 
be  that  initial  value  problems  are  nonself-adjoint  and  thus  they  do  not 
possess  variational  principles  in  the  classical  sense.  In  conjunction 
with  problems  involving  nonconservative  forces,  certain  constrained 
variational  principles  (sometimes  called  extended  Hamilton!s  principles 
—See,  for  example,  ref.  [1])  were  used  for  finite  element  solution 
formulations  [2,  3].  Shortly  afterwards,  using  the  combined  notion  of 
the  Lagrange  multipliers  and  the  adjoint  variable,  some  unconstrained 
variational  statements  were  established  and  used  as  bases  for  finite 
element  solutions  [4,  5].  This  approach  has  been  shown  to  be  more 
advantageous  in  terms  of  simplicity, versatility  and  the  rate  of  conver- 
gence compared  with  the  constrained  variational  approach  [5,  6]. 


1.  Levinson,  M.,  "Application  of  the  Galerkin  and  Ritz  Methods  to  Non- 
conservative Problems  of  Elastic  Stability,"  Zeitschrift  fur 
Angewandte  Mathematic  und  Physik,  Vol.  17,  pp.  431-442  (1966). 

2.  Barsoum,  R.  S.,  "Finite  Element  Method  Applied  to  the  Problem  of 
Stability  of  a Nonconservative  System,"  International  Journal  for 
Numerical  Methods  in  Engineering,  Vol.  3,  pp.  63-87  (1971). 

3.  Mote,  C.  D. , "Nonconservative  Stability  by  Finite  Element,"  Journal 
of  the  Engineering  Mechanics  Division,  Proceedings  of  the  American 
Society  of  Civil  Engineers,  Vol.  EM3,  pp.  645-656  (June  1971). 

4.  Wu,  J.  J.,  "Column  Instability  under  Nonconservative  Forces,  with 
Internal  and  External  Damping  - Finite  Element  Using  Adjoint  Varia- 
tional Principles,"  Development  in  Mechanics,  Vol.  7,  pp.  501-514 
(1973). 

5.  Wu,  J.  J.,  "A  Unified  Finite  Element  Approach  to  Column  Stability 
Problems,"  Development  in  Mechanics,  Vol.  8,  pp.  279-294  (1975). 

6.  Wu,  J.  J.,  "On  the  Numerical  Convergence  of  Matrix  Eigenvalue  Problems 
Due  to  Constraint  Conditions,"  Journal  of  Sound  and  Vibrations, 

Vol.  37,  pp.  349-358  (1974). 
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Fried  was  first  to  treat  the  time-dimension  identically  with  the 
space  dimensions  in  using  the  finite  elements  [7].  His  solution 
formulations,  however,  emanate  from  constrained  variational  principles. 

In  contrast,  this  paper  presents  a generalization  of  the  unconstrained 
variational  approach  to  time-dependent  problems. 

At  this  point,  the  variational  principles  of  integrals  of  convo- 
lution developed  by  Gurtin  [8,  9]  should  be  mentioned.  The  applications 
of  these  principles  in  conjunction  with  finite  elements  in  the  time- 
dimension  [10,  11,  12,  13]  have  so  far  failed  to  show  any  advantage  over 
the  procedure  described  by  Fried.  In  fact,  all  these  analyses  had  to 
resort  to  either  the  Fried's  or  some  other  similar  step-by-step  procedure 
to  complete  the  solutions  in  the  time-dimension. 


7.  Fried,  I.,  '’Finite  Element  Analysis  of  Time  Dependent  Phenomena, ft 
Journal  of  the  American  Institute  of  Aeronautics  and  Astronautics, 

Vol . 7,  No.  6,  pp.  1170-1172  (June  1970). 

8.  Gurtin,  M.  E.  , ,TVariational  Principles  for  Linear  Elastodynamics 
Archive  for  Rational  Mechanics  and  Analysis,  Vol.  16,  No.  1,  pp.  36-50 
(1964). 

9.  Gurtin,  M.  E.,  ”Var5 ational  Principles  for  Linear  Initial-Value 
Problems, " Quarterly  of  Applied  Mathematics,  Vol.  22,  pp.  252-256 
(1964). 

10.  Wilson,  E.  L.  and  Nickel,  R.  E.,  ’’Application  of  the  Finite  Element 
Method  to  Heat  Conduction  Analysis,”  Nuclear  Engineering  and  Design, 
Vol.  4,  pp.  276-286  (1966). 

11.  Dunham,  R.  S. , Nickel,  R.  E.  and  Strickler,  D.  C.,  ’’Integration 
Operators  for  Transient  Structural  Response,”  Computers  and  Struc- 
tures, Vol.  2,  pp.  1-15  (1972). 

12. ^  Ghaboussi,  J.  and  Wilson,  E.  L. , ” Variational  Formulation  of  Dynam- 

ics of  Fluid-Saturated  Porous  Elastic  Solids,:  Journal  of  the 
Engineering  Mechanics  Division,  Proceedings  of  the  American  Society 
of  Civil  Engineers,  Vol.  EM4,  pp.  947-963  (1972). 

13.  Atluri,  S.,  ”An  Assumed  Stress  Hybrid  Finite  Element  Model  for  Linear 
El astodynamic  Analysis,”  Journal  of  the  American  Institute  of  Aero- 
nautics and  Astronautics,  Vol.  11,  No.  7,  pp.  1028-1031  (1973). 
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In  this  paper,  the  use  of  unconstrained  variational  principles  - 
finite  elements  for  usual  boundary  value  problems  is  first  illustrated 
and  the  advantages  over  the  constrained  formulations  are  pointed  out. 

The  unconstrained  variational  principles  can  always  be  constructed 
through  the  use  of  the  Lagrange  multipliers.  The  unconstrained  vari- 
ations are  then  shown  to  lead  naturally  to  (nonself-)  adjoint  variational 
statements.  Thus,  nonconservative  problems  can  be  formulated  easily 
using  finite  elements.  The  application  to  a control  problem  is  given 
[14].  With  the  introduction  of  a cross-product  term  involving  two-point 
boundary  (initial)  values,  the  unconstrained  variational  - finite  element 
formulation  is  again  easily  extended  to  include  time-dependent  problems. 
This  formulation  is  obviously  simpler  compared  with  those  derived  from 
Curtin’s  variational  principles  because  no  convolutional  integrals  are 
needed.  It  is  also  easier  to  use  and  more  versatile  than  the  Fried’s 
procedure  due  to  the  fact  that  no  boundary  or  initial  conditions  are 
involved  in  the  solution  formulation  and  because  of  the  nature  of  the 
Lagrange  multipliers.  As  further  examples  of  application,  finite 
element  matrix  equations  are  derived  for  several  transient  problems 
including  a force  vibration,  a heat  transfer  and  a wave  propagation 
problem.  Detailed  formulations  and  numerical  results  of  two  examples 
are  given  and  comparisons  with  some  known  exact  solutions  are  made. 

2.  LAGRANGE  MULTIPLIERS  AND  FINITE  ELEMENT  FORMULATIONS,  One  of 
the  advantages  of  the  finite  element  method  is  its  capability  of  solving 
large  complicated  problems  in  a routine  manner.  However,  the  same  con- 
cepts used  in  a program  for  large  systems  may  be  understood  using 
relatively  simple  problems. 

Let  us  consider  the  stability  of  a Euler’s  column.  The  governing 
equations  are  as  follows: 


14.  Wu,  J.  J.,  "On  the  Stability  of  a Free-Free  Beam  under  Axial  Thrust 
Subjected  to  Directional  Control,"  Journal  of  Sound  and  Vibration, 
Vol.  43,  pp.  45-52  (1975).  Also  see  a correction  note  on  this  paper, 
ibid,  Vol.  44,  p.  309  (1976). 


3 


D.E. 

E I u""  + P u"  + co2  pAn  = 0 

(la) 

B.C. 

u (0)  = u' (0)  = 0 

(lb) , (lc) 

o 

It 

t — 

o* 
v / 

'h 

• (Id) 

E I u"'  (Z)  +Pu'(l)  = 0 

(le) 

A usual  variational  principle  can  be  written 

6Jj  (u)  = 0 

(2a) 

where 

1 n 

JlCu) 

= 1 4 Cu")2  - P(u')2  + u2pAu2]  dx 

(2b) 

To  establish 

the  equivalence  between  eqs.  (1)  and  (2), 

one  simply 

carries  out  the  variation  of  Jj  in  eq.  (2a): 

6Ji  = 

Jq  [E  I u"6u"  - P u'6u'  + co2  pA  u <5u]  dx 

(3a) 

= Iq  [E  I u""  + P u"  + a)2  pAu]  6u  dx 

+ 

[E  I u"  6u'  -(El  u'"  + P u')  6u]x  _ o 

- 

[E  I u"  6u'  - (El  u’"  + P u')  6u]x  = 0 

(3b) 

From  eq.  (3b)  one 

observes  that  for  the  coordinate  functions 

and  their 

variations  satisfying  the  boundary  conditions  in  eqs.  (lb  - le) , eq. 
(la)  implies  eq.  (2a)  and  vice  versa.  The  finite  element  formulation 
for  this  problem  begins  with  eq.  (3a). 


Let 


u(x)  = a1  (x)  U 


(4) 


where  a(x)  is  the  displacement -function  vector  and  U , the  generalized 
displacement  vector.  Upon  the  substitution  of  eq.  (4)  into  eq.  (3a), 
one  immediately  obtains 


5UJ 


| Kx  + w2M  J.  U = 0 


where 


Kj  = J*  [ E I a"  a"T  - Pa'  a,T]  dx 


,.Ti 


M = f pA  a aT  dx 


(5) 

(6a) 

(6b) 
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Eq.  (5)  is  not  yet  ready  to  be  solved  since  neither  U nor  6U 
consists  of  independent  elements  due  to  the  boundary  conditions  re- 
quirements placed  on  u(x). 

Let  us  now  consider  a slightly  different  variational  principle: 

6J2  = 0 (7a) 

with 

J2  = ~ [ El  (u")2  - P (u')2  + O)2  pAu2]  dx 

+ I ax  [u(0)] 2 + j a2  [u*  CO)] 2 (7b) 

where  and  are  the  Lagrange  multipliers. 

Carrying  out  the  variation  of  eqs.  (7),  we  have 
£ 

6J~  = f [El  u"  6u"  - Pu'  6u  + w2pAu2]  dx 
1 o 

+ u(0)  6u(0)  + a 2 u'  (0)  Su'(0)  (8a) 

= y [E  I uM"  + Pu"  + a)2  pAu]  6u  dx 
o 

+ [E  I u"  <Su'  - (E  I u'"  + Pu')  6u]  * 

- [(E  I u"  - a^u'  -(El  u1"  + Pu'  + a]u)  6u]  ^ (8b) 

Eq.  (8b)  states  that  a necessary  and  sufficient  condition  for 
6J2  = 0 is  the  problem  defined  by  the  following  set  of  equations: 


E I n""  + 

Pu"  + to2pAu  = 0 

(9a) 

E I u" (0) 

- a2  u'  (0)  = 0 

(9b) 

C, 

/ — \ 
O 

V — / 

•f  ■ 

Pu'  (0)  + ax  u(0)  = 0 

(9c) 

E I 

C 

v — ' 
ii 

O 

(9d) 

E I u"' 

(2,)  + Pu'  (2.)  = 0 

(9e) 

provided  that  the  variation  Su  is  completely  arbitrary,  comparing 
eqs.  (9)  and  (1),  it  is  seen  that  eqs.  (1)  is  a special  case  of  (9)  as 
, a2  approach  to  infinity.  From  eq.  (8a),  we  can  see  that  the  finite 
element  matrix  equation  now  becomes 
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where 


(10) 


(11) 


The  matrix  JC  in  eq.  (11)  has  been  defined  in  eq.  (5)  and  the  super- 
script T denotes  the  transpose  of  a matrix  (a  vector) . Since  6u  is 
arbitrary,  6\J  in  eq.  (10)  is  arbitrary,  eq.  (10)  leads  directly  to  the 
final  matrix  equation  to  be  solved. 


It  is  then  clear  that  the  method  of  Lagrange  multipliers,  used  in 
conjunction  with  the  finite  element  method,  will  not  only  facilitate 
the  solution  formulations  but  also  encompass  a larger  class  of  problems 
to  be  solved  compared  with  the  use  of  constrained  variational  statements. 
The  applications  of  the  same  general  concept  can  be  extended  further. 

3.  FROM  UNCONSTRAINED  VARIATIONS  TO  ADJOINT  VARIATIONAL  STATEMENTS. 

We  have  noted  that  the  variation  6u  in  eq.  (8)  is  quite  independent 
of  the  function  u itself  and  nothing  will  be  changed  if  we  simply 
replace  <5u  with  6v  to  emphasize  this  independence.  This  substitution, 
however,  has  suggested  the  adjoint  variational  principles.  Let  us 
consider 


(12) 


6J3  = 0 
Si 

Jj  = f (El  u"v"  - Pu'v'  + u2pAuv)  dx 
+ a1u(0)v(0)  + a2u'(0)v(0)  + a3Pu'  (£)v(£) 


(13b) 


(13a) 


Carrying  out  the  variations,  we  have: 


«J3  - (6J3)u  . («J3)v 


(14) 


where 


u o 


+ u(0)<5v(0)  + o^u1  (0)6v' (0)  + a^u*  (X)  Sv(£) 


(15a) 


+ [E  I u"6v'  -(El  u'"  + Pu'  - a3u')  6v] 

X — )C 
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- [ (E  I u"  - a2u')  6v'  - (E  I u'"  + Pu'  + cijU)  6v]x  = Q (15b) 
and 

(6J3)v  = J*  (E  I v"  <5u"  - Pv'6u'  + (jjzpAv6u)  dx 

+ ai  v (0)  6u(0)  + a2  v'  (0)  5u'  (0)  + a3  v(Jl)  6u'  (£)  (16a) 

l 

= J (El  v""  + Pv"  + u2pAv)  6u  dx 
o 

+ [(E  I v"  + a3v)  6u'  -(El  v"'  + P v')  6u] 

X " jC 

- [ (E  I v"  - a?  v'.)  6u'  - (El  v"'  + Pv'  + ai  v)  6ul  (16b) 

* 1 x = 0 

From  eq.  (15a),  it  is  clear  that  a necessary  and  sufficient 
condition  for  (SJ3)U  =0  is  the  problem  defined  by  the  following  set 
of  equations: 


D.E.  E I u,m  + Pu"  + oo2pAu  = 0 (17a) 

B.C . E I u"(«,)  = 0 (17b) 


E I u'"  («,)  +(P  - a3)  u'  (£)  = 0 


(17c) 


E I u" (0)  - a2  u'  (0)  = 0 


(17d) 


E I u'"  (0)  + P u'(0)  + ax  u(0)  = 0 (17e) 

Now  eqs.(9)  has  become  a special  case  of  eqs.(17)  when  a3  = 0. 

In  addition,  the  problem  defined  by  (<$J3)V  = 0 of  eqs.  (16)  is  called 
the  adjoint  problem  to  eqs.  (17).-  For  a3  = 0,  the  adjoint  problem  is 
identical  to  the  problem  itself  — hence,  the  self-adjoint  system.  Now, 
considering 

a3^=  k P (18) 


in  eq.  (17c) , we  have 


E I u'"  (£)  -K  P u'  (£)  = 0 


(19) 


K = k - 1 


(20) 
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Eq.  (19)  defines  the  boundary  condition  of  a general  non-conser- 
vative load.  It  is  also  clear  from  eq.  (19)  that  K is  a dimensionless 
design  constant  which  defines  the  small  angle  between  the  direction  of 
the  applied  load  P and  the  tangent  of  the  deflected  column  at  the  end. 
Since  (SJ3)U  = 0 alone  defines  the  boundary  value  problem  of  eq.  (17) 
and  vice  versa,  we  need  not  at  all  to  be  concerned  with  the  adjoint 
problem.  Now  it  is  a simple  matter  to  modify  the  finite  element  matrix 
equation  as 

| K3  + u2M  | U = ° (22) 

where 

K3  = K2  + Ct3  a’  (£)  aT(£)  (23) 

4.  FINITE  ELEMENTS  FOR  INITIAL  AND  INITIAL- BOUNDARY  VALUE  PROBLEMS. 

(1)  A Forced  Vibration  Problem.  Let  us  first  consider  a problem 
of  "one"  degree  of  freedom,  i.e.,  a mass-spring  system.  The  differential 
equation  and  initial  conditions  are 


m u + k u = f(t),  0<t<T 

(24a) 

u (0)  = uQ 

(24b) 

u(0)  -Uj 

(24c) 

where  u(t)  is  the  displacement  of  the  mass  centre  from  its  equilib- 

rium position,  m , the  amount  of  mass  and  k , the  spring  constant. 
The  function  f(t)  is  given,  so  are  the  constants  a and  b.  The  constant 
T appeared  in  the  bounds  of  eq.  (24a)  is  any  given  positive  number 
other  than  infinity.  In  order  to  formulate  approximate  solutions  for 
eqs.  (24)  the  way  we  did  in  the  previous  section,  let  us  consider  a more 
general  case 

m ii  + k u = f (t)  ^ 

(25a) 

u(T)  - a [ u(0)  - uQ  ] =0 

(25b) 

u(0)  = u2 

(25c) 

where  a is  a parameter,  obviously  eqs.  (25)  reduce  to  (24) 
a approaches  to  <»  . Now,  with  eqs.  (25),  we  are  able  to  write 
unconstrained  variational  statement  as  follows: 

when 

an 

6 J4  = 0 

(26a) 
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where 


T 

JA  = f [ - muv  + kuv  - f(t)v  ] dt 
4 o 

+ ma  [u(0)  - Uq]  v(T)  - mu^  v(0)  * (26b) 


Since 

T 

($J4)  = J [ - mufiv  + ku6v  - f(t)  6v  ] dt 

4 U Q 

+ ma  [ u(0)  - Uq  ] 6v(T)  - mu^6v(0)  (27a) 

= jT  [ mu  + ku  - f(t)  ] <5v  dt 

- m | u(T)  - a [u(0)  - uQ]  j <5v(T) 

+ m [u(0)  - uj  ] 6v(0)  (27b) 

The  already  familiar  form  of  eqs.  (27)  state  that  (a),  (<SJ)U  = 0 is  a 
necessary  and  sufficient  condition  for  eqs.  (25), and  (b),eq.  (27a) 
provides  us  the  finite  element  matrix  equation.  Thus,  if  we  assume  as 
before  that 


u(t)  = aT(t)  U 

v (t)  = aT(t)  V 

t 

Eq.  (27a) 

yields 

6VT  k4  U = 6VT  f 

(28) 

where 

T 

K4  = /0  ( -m  a aT  + k a aT)  dt 

+ ma  a(T)  a^(0) 

(29) 

and 

F 

T 

= J f(t)  a dt  + mauQ  a(t)  + m uq  a(0) 

(30) 

Again,  since  6V  is  unconstrained  eq.  (28)  leads  directly  to 

K4  U = F (31) 

which  is  the  final  equation  to  be  solved. 
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(2)  A Heat  Conduction  Problem.  The  one  dimensional  transient 
heat  conduct  problem  can  be  described  by  the  equation 


3 3u  3u 

3x  (K  3x  pc  3t  “ = 0 

Boundary  and  initial  conditions  are 

u(0,t)  = g0(t) 

= g (t) 

u(x,0)  = h (x) 

K = thermal  conductivity 

p = material  density 

c = specific  heat 

f(x,t)  = heat  source  function 

gQ(t),  glOO  anc^  h(x)  are  prescribed  functions 

Let  us  consider 

<$J5  - 0 


where 


and 


since 


JS  = -J,LJT  v]  dt  dx 

+ /o  aK  [u(L,t)  - gx Ct) ] v(L,t)  dt 

- JoT  aK  [u(0,t)  - g0(t)]  v (0,t)  dt_ 

,L 

- J pc  [u(x,0).  - h(x)]  v(x,0)  dx 

0 

(<5J5)U  “jV*  P-  + pc  fj.  6v  + f(x,t)  6v]  dx  dt 

o'o  oX  oX  ot 

+ /oT  a K [U (L,t)  - gl(t)]  <5v(L,t)  dt 

/0T  a K [u(0,t)  - g0(t)]  6v(0,t)  dt 
, L 

J pc  [u(x}0)  - h(x)]  5v(x,0)  dx 


(32a) 

(32b) 

(32c) 

(32d) 


(33a) 


(33b) 


(34a) 
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- i'll  57  (K  s ) - PC  !r  - f(x,t)]«v  dxdt 


f 


+ 


3u(L,t) 
1 9x 

8u(0,t) 

3x 


L 


- a[u(0,t)-g1 
a[u(0,t)  - gQ 


+ Jo  pc  [u(x,0)  - h(x)  ] 


(t)]  J-Sv(L,t)  dt 

(t)]  |>6v (0,t)  dt 

6v(x,0)  dx 


(34b) 


it  is  clear  that  (6J5)U  = q is  a necessary  and  sufficient  condition 
for  eqs  (32)  as  a °o  ana  eq.  (34a)  provides  the  finite  element  ma- 
trix equation.  We  can  write  from  eq.  (34a), 


- fr.V&'Hg)  * P'  H !'l  dxdt 

,T 

+ jo  aK  [u(L,t)  Sv(L,t)  - u(0,t)  <5v(0,t)]  dt 
L 

+ j p c u(x,0)  6v(x,0)  dt 
L T 

= Jo  /0f  (x>t)  6v  dxdt 

T 

+ Jo  aK  [g1(t)  6v(L,t)  gp(t)  Sv(0,t)]  dt 

L 


+ J pc  h(x)  6v(x,0)  dx 

(35) 

Now,  let 

u(x,t)  = aT(x,t)  U 

(36a) 

v(x,t)  = a^(x,t)  V 

(36b) 

in  the  usual  manner,  we  have 

SVT  K U = 6VT  F 

(37) 

Li  T 

K = -*  f f (K  a Y at,  + p c a at)  dx  dt 

~ JoJo  ~ > X ~ ~ > L 

T 

+ / aK  [a(L,t)  a(L,t)  - a(0,t)  aT(0,t)]  dt 


L 

+ / p c a(x,0)  a^(x,0)  dx  (38) 

0 • ~ - 
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and 


LT 

F = U0fCx,t)  a(x,t)  dxdt 
T 

+ fQ  olK  [gl(t)  a(L,t)  - gQ(t)  a(0,t)]  dt 
,L 

+ J0  p c h(x)  a(x,0)  dx  (39) 

Again,  since  SV  in  eq.  (37)  is  completely  arbitrary,  we  arrive  at  the  •• 

final  matrix  equation  to  be  solved. 

K U = F (40) 


(3)  A Wave  Propagation  Problem.  For  a quite  general  wave  propa- 
gation problem,  the  following  system  can  be  written. 


9ni 

9xJ 


- c 


2 9fu 
dt2 


= f(x,t) 


(41a) 


u(0,t)  = gn(t) 


u(L,t)  = gl(t) 


u(x,0)  = hQ(x) 
u(x,0)  = h:(x) 


(41b) 
(41c) 
(4  Id) 
(41e) 


The  extension  of  the  previous  formulation  to  this  problem  is  straight 
forward.  Let  us  consider 


where 


-L.T 


«J6  • 0 


t - r 9u  9v  . 2 9u  9v  ct 

J6"U0  [-9T9x  + C 9?  9t  £Cx>t)  V]  dxdt 


(42a) 


- a /0  [u(L,t)  - g1(t)]  v (L,t)  dt 

J 

+ a /q  [u(0,t)  - gQ(t)]  v (0,t)  dt 
L 

- a Jo  [u(x,0)  - hQ(x)]  v (x , T)  dx 


+ J0  [u (x, 0)  - h1(x)]  v(x,0)  dx 


(42b) 
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Again, 


C<SJ6)U  = /oL/oT  [ 6(  |r  5 + c2  §£  *(  §£  ) - f(x,t)6v]  dxdt 


T 

-a  / [u (L,t)  - g.(t)]  6v(L,t)  dt 

Q 1 

,t 

+ a Jo  [u (0,t)  - gQ(t)]  5v(0,t)  dt 
rL 

- a Jo  [u(x,0)  - hQ(x)]  6 v(x,T)  dx 


+ J [u 
0 

(x,0) 

- h1(x)]  6v(x,0) 

dx 

(43a) 

LT  ' 32, 

:U  t 8 ’ 

u 

C 

2 3zu 

f(x,t)]  V 

dxdt 

0 0 3x 

2 

9t2 

+ 4 i 

f 3u(L,t) 
[ 3x 

- a 

[u (L,t) 

- gi(t)]  } 

6v(L,t) 

dt 

-4 1 

f 3u(0,t) 

[ 3x 

- a 

[u(0,t) 

- SqM]  } 

Sv(0,t) 

dt 

Hi 

r 3u(x,t) 

[ 3t 

- a 

[u (x,0) 

- h0W3  } 

6v(x,T) 

dx 

-/  t 

V 

3u , 
— (x 

,0)  - hx 

(x)]  <5v(x,0)  dx 

(43b) 

From  eqs.  (43),  it  is  again  clear  that  (6Jg)u  =0  is  a necessary  and 
sufficient  condition  for  eqs.  (41)  as  a -*■  00  and  that  eq.  (43a)  will 
yield  the  finite  element  matrix  equation.  From (43a)  one  has: 

L T 

6(|1)+  c2  |H.6(|1)  \ dx  dt 
J0  Jo\  dx  KdxJ  3t  '•at*'  J 

- a J u(L,t)  fiv(L,t)  dt  + a J u(o,t)  6v(0,t)  dt 

- a f u(x,0)  6v(x,T)  dx 

v 

= ff  f(x,t)  5v(x,t)  dx  dt 
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Again,  let 


(44) 


a Jt  6 v (L,t)  dt  + a 4 g0ft)  6v(0,t)  dt 

- a 4 Vx)  6v(x,T)  dx  + 4 hlW  6v (x,0)  dx 


u(x,t)  = a (x,t)  U 

(45a) 

v(x,t)  = aT(x,t)  V 

(45b) 

Eq.  (44)  becomes. 

in  matrix  form. 

6 VT  K U = 6 VT  F 

(46) 

where 

T L ^ 

$ ” 44  ft  2*  2*^  + a ah  dxdt 

T rp  T T 

■ a 4 sftft)  t ft,t:) dt  + a 4 aT(°ft) dt 


- a 4 a(x,t)  aT(x,0)  dt  (47) 

T L 

P * 44  f(xft)  aftft)  dxdt 

T T 
- a j0  gx(t)  a(L,t)  dt  + a /o  gQ(t)  a(0,t)  dt 

-L  L 

+ “ 4 h0°°  2(XjT)  dt  + 4 bjW  a(x,°)  dt  (48) 

Due  to  the  arbitrariness  of  6V  , eq.  (46)  leads  directly  to  the  final 
matrix  equation 


£H=£  (49) 

5,  NUMERICAL  DEMONSTRATIONS,  Several  numerical  examples  will  be 
given  in  this  section  to  demonstrate  the  application  of  the  formulation 
described  so  far. 

(!)  Forced  Vibration.  We  shall  consider  a special  case  of  the 
forced  vibration  problem  formulated  earlier.  The  forcing  function  in 
eqs . (24)  is  taken  to  be  a cosine  function  thus,  rewrite  eqs.  (24), 

mu+ku=f0  cos  wft  (50a) 
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u(0)  = u0  (50b) 

u(0)  = u:  (50c) 

where  Uq,  u^,  fp  and  co£  are  given  constants.  In  the  finite  element 
formulation,  we  shall  replace  eqs.  (50)  with  the  following  set 

m il  + k u = fQ  cos  u)£t  (51a) 

u(t)  - a [u(0)  - u0]  = 0 (51b) 

u(0)  - Uj  = 0 


thus,  eqs.  (50)  becomes  a special  case  of  (51)  as  a 00 . It  is  con- 

venient to  nondimensionalize  the  independent  variable  t and  let 


x = t/T 


In  terms  of  x,  eqs.  (51)  become 


where 


..  2 ? 

u + T 0)  u = cos  (T  (i3f  t) 
u(l)  - T a [u(0)  - Uq]  = 0 
u(0)  - T Ul  = 0 
fj  = T2  fQ/m  oo2  = k/m 


The  exact  solution  for  eqs.  (53)  can  be  easily  written  as 
u(x)  = A cos  (T  w • x)  + B sin  (T  to  • x) 
+ q cos  (T  Wf  • x) 


with 


Ui 


n = 


— - — , e - 

m((o2-Wf2)  w 


(52) 

(53a) 

(53b) 

(53c) 

(54) 


(55) 


a u0  t T uj  cos  (Too)  - r\  [a  + Tto£  sin  (Tco£)] 
a + T to  sin  (T  co) 


(56) 


To  solve  eqs.  (53)  using  finite  elements,  one  begins  with  the  variational 
statement: 


6 J » 0 

J = / [-  u v + T2  (o2uv  - f (x)  v]  dx 


(57a) 
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Now  that 


(57b) 


+ T a [u(0)  - u0]  v(l)  -Tuj  v(0) 


(6J)U  = 0 


From  eq. 


with 


eq.  (59) 
or 

where 


1 

= Jo  [-  u 6 v + T2  a)2  uSv  - f(x)6  v] 

+ T a [u(0)  - uQ]  v(l)  - TUJ6VC0) 

= Jo  [u  + T2  co2  u - £ (t)  ] <5v  dt 

- (u(0)  - a T[u(0)  - u0]}  <5v(1) 

+ { u(0)  - T Ul>  6v(0) 

(58b) , one  has  . 

1 

J0  [-  u 6 v + T2w2  u 6 v]  dt  + aT  u(0)  6v(l) 

1 

“ J0  f(T)  6v  dx  + a T uQ  <$v(l)  + T Uj  6v(0) 
u(x)  = aT(x)  U 
v(x)  = aT(x)  V 

leads  to 

6 VT  K U = 6 VT  F 
K U = F 

K = / (-  4 aT  + T2  u2  a aT)  dt 
+ a T a(l)  a (0) 

1 

F = / f(x)  a dx  + aTu0  a(l)  + T u:  a(0) 


(58a) 


(58b) 


(58c) 


(59) 

(60) 


(61) 

(62) 

(63) 


The  results  obtained  from  this  finite  element  formulation  are  com- 
pared with  the  exact  solutions  as  shown  in  Tables  1-3.  The  values  of 
the  parameters  chosen  for  these  data  are  k=1.0,  m=1.0,  fg  * 1.0, 

0)£  = 0.5,  uq  = 1.0,  Uq  = 1.0  the  number  of  elements  used  is  ten.  The 
calculated  u and  u for  T = 2.0,  10.0  and  20.0  are  given  in  Table  1, 

2,  3 and  4 respectively.  The  forcing  function  cos  a)ft  and  the  solu- 
tion u(t)  are  also  plotted  in  the  range  0 < t < 20  as  shown  in  Figure  1. 

(2)  Solutions  to  a Transient  Heat  Conduction  Problem.  As  another 
numerical  example,  we  shall  take  the  nondimensional  heat  transfer  problem 
defined  by  the  following  set: 
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TABLE  1.  Solutions  to  the  Forced  Vibration  Problem  Using 

FE-UVF  Compared  with  Exact  Solutions  (in  Parentheses) 

0 < t < 2.0 

^ c 
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TABLE  2.  Solutions  to  the  Forced  Vibration  Problem  Using 


FE- 

UVF  Compared  with 
0 $ t * 

Exact  Solutions 
10.0 

(in  Parentheses) 

t 

u(t) 

u(t) 

0 

1.000 

(1.000) 

1.004 

(1.000) 

1.0 

1.832 

(1.831) 

0.505 

(0.501) 

2.0 

1.770 

(1.768) 

- 0.675 

(-0.674) 

3.0 

0.566 

(0.565 

- 1.614 

(-1.608) 

4.0 

- 1.094 

(-1.094) 

- 1.518 

(-1.512) 

5.0 

- 2.123 

(-2.122) 

- 0.435 

(-0.435) 

6.0 

- 1.920 

(-1.919) 

0.778 

(0.773) 

7.0 

- 0.843 

(-0.843) 

1.213 

(1.207) 

8.0 

0.167 

(0.166) 

0.690 

(0.689) 

9.0 

0.436 

(0.435) 

- 0.126 

(-0.122) 

10.0 

0.114 

(0.114) 

- 0.385 

(-0.381) 
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TABLE  3.  Solution  to  the  Forced  Vibration  Problem  Using 

FE-UVF  Compared  with  Exact  Solutions  (in  Parentheses) 
0 < t £ 20.0 


t 

u(t) 

u(t) 

- 0 

1.000 

(1.000) 

1.05 

(1.00) 

2.0 

1.778 

(1.768) 

- 0.68 

(-0.67) 

4.0 

- 1.097 

(-1.094) 

- 1.57 

(-1.51) 

6.0 

- 1.928 

(-1.919) 

0.82 

(0.77) 

8.0 

0.173 

(0.166) 

0.71 

(0.69) 

10.0 

0.116 

(0.114) 

- 0.44 

(-0.38) 

12.0 

0.453 

(0.462) 

0.88 

(0.85) 

14.0 

1.956 

(1.950) 

0.06 

(0.03). 

16.0 

- 0.156 

(-0.162) 

- 1.76 

(-1.71) 

18.0 

- 2.199 

(-2.186) 

0.15 

(0.14) 

20.0 

- 0.348 

(-0.342) 

1.10 

(1.08) 
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Figure  1 . Forcing  Function  F(t)  and  Solution  u(t) 
for  the  Vibration  Problem 
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S^u  £u 

D.E. : 

r “ = 0 , 0<x<l  ; 0<t<T 

3x2  3t 

(64a) 

B.C . : 

u(0),t)  - 1,  |!L 
3x 

(l,t)  = 0 

(64b, c) 

I.C.: 

u(x,0)  = 

0 

' (64d) 

where  T is  any  given  finite  real  positive  number.  To  facilitate  compu- 
tation, it  is  desirable  to  change  the  independent  variable  t into  T 
such  that 

T = t/T  (65) 


thus,  the  system  of  eqs.  (64)  becomes 


D.E.  : 

— t - — = 0 , 0<x<l  ; 0<t<1 
3x2  T 3t 

(66a) 

B.C.  : 

u (0 , t)  = 1 ; (1 , t)  = 0 

3x 

(66b, c) 

I.C.: 

u(x,0)  = 0 

(66d) 

According  to  our  unconstrained  variational  formulation, 
again  replaced  by  the  following: 

this  system  is 

D.E.  : 

2 - T = 0 , 0<x<l  ; 0<t<1 

3x2  " aT 

(67a) 

B.C. : 

3u 

(0,t)  + a [u(0,t)  -1]  = 0 

(67b) 

§■  d.T)  = 0 

(67c) 

I.C. : 

u(x,0)  = 0 

(67d) 

Clearly,  eqs.  (67)  reduces  to  (66)  as  a 00  . The  variational  state- 
ment can  be  written  as 

6J  = 0 (68a) 

where  1 

J - - / 

'o 

1 

+ ofo  [ll(0,T)  - 1]  V(0,T)  dT 
1 

+ /0  u(x,0)  v(x,0)  dx  (68b) 

Due  to  the  fact  that  v(x,t)  is  unconstrained,  it  is  a simple  matter  to 
show  that 

(6J)U  « 0 (69) 


J1  ( 9u  3v  + I 9u  y)  dxdt 
o 3x  3x  T 3x 
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is  a necessary  and  sufficient  condition  of  eqs . (67)  . Now  the  finite 
element  matrix  equations  can  be  obtained  from  eq.  (69). 

(W)u  = " /V  [—  «(— ) + - — 6v]  dxdt 

u 0 'o  3x  3x  T 9t 

1 

+ / [u(0,x)  - 1]  6v(0,x)  dt 

1 

+ J u(x,0)  6v(x,0)  dx  = 0 (70) 

o 

or, 

- jV  [&L  6(~)  + 1 2"  6v]  dxdt 

o o 3x  3x  9x 

, i 1 

4 a f u(0,x)  6v(0,t)  dt  4 f u(x,0)-  6v(x,0)  dx 
o o 

1 

= a J 6v(0,t)  dt  (71) 

Using  the  usual  procedure  of  discretization  and  the  assumption 
of  displacement  functions,  the  final  finite  element  matrix  equation 
evidently  can  be  derived  from  eq.  (71).  We  shall  omit  the  details  here. 
The  computational  results  are  presented  in  Tables  4 and  5.  The  finite 
element  grid  scheme  used  is  shown  in  Figure  2.  As  clearly  shown  in 
those  tables,  excellent  agreement  exists  between  the  FE-UVF  approach 
and  the  series  solution.  It  is  noted  that  the  approximate  solutions 
are  less  accurate  invariably  as  they  approach  the  initial  time  t = 0. 
This  is  probably  due  to  the  discontinuity  of  the  initial  boundary  data 
at  x = 0,  t = 0. 
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TABLE  4.  Transient  Heat  Transfer  Solutions  u(x,t)  Using  FE-UVF 
Compared  with  Exact  Series  Solutions  (in  Parentheses) 

0<t<T=1.00 


V X 

t N, 

0 

0.2 

0.4 

0.6 

0.8 

1.0 

0.2 

1.000 

0.754 

— 

0.583 

0.370 

0.264 

0.228 

(1.000) 

(0.757) 

(0.496) 

(0.405) 

(0.284) 

(0.179) 

0.4 

1.000 

0.855 

0.713 

0.622 

0.552 

0.516 

(1.000) 

(0.853) 

(0.721) 

(0.616) 

(0.549) 

(0.526 

0.6 

1.000 

0.910 

0.828 

0.767 

0.725 

0.708 

(1.000) 

(0.910) 

(0.830) 

(0.767) 

(0.724) 

(0.710) 

0.8 

1.000 

0.945 

0.896 

0.857 

0.832 

0.823 

(1.000) 

(0.945) 

(0.896) 

(0.857) 

(0.832) 

(0.823) 

1.0 

1.000 

0.967 

0.937 

0.913 

0.897 

0.892 

(1.000) 

(0.967) 

(0.937) 

(0.913) 

(0.897) 

(0.892) 
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TABLE  5.  Transient  Heat  Transfer  Solutions  u(x,t)  Using  FE-UVF 
Compared  with  Exact  Series  Solutions  (in  Parentheses) 

0<t<T=0.05 


X 

t . 

0 

0.2 

0.4 

0.6 

0.8 

1.0 

0.01 

1.000 

0.144 

0.014 

0.002 

0.000 

0.000 
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Figure  2.  Finite  Element  Grid  Scheme  Used  for  a 
Transient  Heat  Conduction  Problem 
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